Code
cor(pooled_data$p_vra, pooled_data$ssq_full)[1] 0.3877125Executive Summary
This report investigates the association between pre-VR conditions and change in total symptom score before and after undergoing VR. This analysis specifically focuses on VR anxiety measured by p_vra, the expectancy score measured by p_e, and the change measured by ssq_full. The results suggest an association between these variables.
Explorary Data Analysis (EDA)
Structures
Limitations
Assumptions
Research Question and Analysis
Research Question
RQ1 : Is there an association between anxiety before undergoing VR and change in total symptom score?
Analysis
Code
figure_1 <- ggplot(pooled_data, aes(x = p_vra, y = ssq_full)) +
geom_point() +
geom_smooth(method = "lm", se = FALSE) +
labs(title = "ssq_full vs p_vra")
ggplotly(figure_1)`geom_smooth()` using formula = 'y ~ x'Code
model = lm(ssq_full ~ p_vra, data = pooled_data)
figure_2 <- ggplot(model, aes(x = .fitted, y = .resid)) + # Plotting fitted values and residuals
geom_point() + # Create scatterplot
geom_hline(yintercept = 0, linetype = "dashed", colour = "red") +
labs(title = "Residual Plot")
ggplotly(figure_2)Code
summary(model)
Call:
lm(formula = ssq_full ~ p_vra, data = pooled_data)
Residuals:
Min 1Q Median 3Q Max
-32.977 -10.385 -3.590 7.204 86.410
Coefficients:
Estimate Std. Error t value Pr(>|t|)
(Intercept) 3.7997 1.4447 2.630 0.00893 **
p_vra 2.5968 0.3378 7.687 1.7e-13 ***
---
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Residual standard error: 15.75 on 334 degrees of freedom
Multiple R-squared: 0.1503, Adjusted R-squared: 0.1478
F-statistic: 59.09 on 1 and 334 DF, p-value: 1.697e-13The correlation between p_vra (VR anxiety) and ssq_full (Full Sum-Scored Active SSQ minus Full Sum-Scored Baseline SSQ) shows 0.3877125, indicating low positive correlation. Figure 1 supports linear relationship between two variables. Figure 2 reinforces the relationship by showing approximately random and same variances in vertical direction along the fitted axis (homoscedasticity).
With these assumptions, performing regression test is approriate. The regression test in the summary statistics shows p value of 1.7e-13 and test statistic of 7.687. Assuming the significant level is 0.05, the null hypothesis (the slope is 0) can be rejected, suggesting the alternative hypothesis (the slope value is significantly different from 0). Thus, the alternative hypothesis suggests a positive linear relationship between VR anxiety and change in total symptom score.
Code
PSTAI_vars <- c("PSTAI_1", "PSTAI_2", "PSTAI_4", "PSTAI_6")
median_summary <- lapply(PSTAI_vars, function(var) {
pooled_data %>%
mutate(PSTAI_binary = case_when(
.data[[var]] %in% c("Very Much", "Moderately", "Somewhat") ~ "Positive",
.data[[var]] == "Not at all" ~ "Negative",
TRUE ~ NA_character_
)) %>%
group_by(PSTAI_binary) %>%
summarise(median_ssq = median(ssq_full, na.rm = TRUE)) %>%
mutate(Variable = var)
}) %>%
bind_rows()
median_summary <- median_summary %>%
mutate(Label = recode(Variable,
"PSTAI_1" = "Calm",
"PSTAI_4" = "Relaxed",
"PSTAI_2" = "Tension",
"PSTAI_6" = "Worry"))
median_summary <- median_summary %>%
mutate(Label = factor(Label, levels = c("Calm", "Relaxed", "Tension", "Worry")))
ggplot(median_summary, aes(x = Label, y = median_ssq, fill = PSTAI_binary)) +
geom_col(position = position_dodge(width = 0.7), width = 0.6) +
scale_fill_manual(values = c("Positive" = "#377eb8", "Negative" = "#e41a1c")) +
labs(title = "Median 'ssq_full' Scores by PSTAI Group",
x = "PSTAI Variable",
y = "Median 'ssq_full' Score",
fill = "PSTAI Group") +
theme_minimal()
PSTAI data on the 4-point Likert scale were aggregated into two categories: “Positive” (“Very much”, “Somewhat”, “Moderately”) and “Negative” (“Not at all”). A comparative bar chart showed that change in total symptom score aligned with anxiety, as indicated by higher reports of tension and worry and lower reports of relaxation and calmness, consistent with established markers of anxiety (American Psychiatric Association, 2022).
This is evidenced by the median scores of ssq_full, where the ssq_full revealed higher median “Negative” scores for Calm and Relaxed, and higher median “Positive” scores for Tension and Worry.
Code
median_expect <- median(pooled_data$expect, na.rm = TRUE)
pooled_data$expect_group <- ifelse(pooled_data$expect <= median_expect, "Low", "High")
p <- ggplot(pooled_data, aes(x = expect_group, y = ssq_full, fill = expect_group)) +
geom_boxplot() +
labs(
title = "Symptom Score by Expectancy Group",
x = "Expectancy Group",
y = "Symptom Score (ssq_full)"
) +
theme(
plot.title = element_text(hjust = 0.5)
)
interactive_plot <- ggplotly(p)
interactive_plotWhen looking at the box plot above, the High expectancy group median generally shows higher ssq_full outcomes compared to the Low expectancy group median. Expectancy was divided by whether scores were above or below the overall median. One noticeable feature is an outlier in the Low group, which may explain why the difference between the two groups was less clear in the first analysis. The High group also has a wider IQR, pointing to greater variation in symptom severity among those with higher expectancy.
Taken together, these findings suggest that higher expectancy is linked to both stronger symptom reporting in VR and more variability across participants. These may suggest people who have high expectancy or negative expectancy have anxiety (Steinman et al. 2013), resulting in greater change in total symptom score.
Articles
Code Chunks. (n.d.). Rmarkdown.rstudio.com. https://rmarkdown.rstudio.com/lesson-3.html
Title Blocks – Quarto. (2018). Quarto. https://quarto.org/docs/authoring/title-blocks.html
American Psychiatric Association. (2022). Diagnostic and statistical manual of mental disorders (5th ed., text rev.). Arlington, VA: American Psychiatric Publishing.
Steinman, S. A., Smyth, F. L., Bucks, R. S., MacLeod, C., & Teachman, B. A. (2013). Anxiety-linked expectancy bias across the adult lifespan. Cognition & Emotion, 27(2), 345–355. https://doi.org/10.1080/02699931.2012.711743